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We consider nonlinear modes of the nonlinear Schrodinger equation with a nonlocal nonlinearity 
and an additional PT-symmetric parabolic potential. We show that there exists a set of continuous 
families of nonlinear modes and study their linear stability in the limit of small nonlinearity. It 
is demonstrated that either VT symmetry or the nonlocality can be used to manage the stability 
of the small-amplitude nonlinear modes. The stability properties are also found to depend on the 
particular shape of the nonlocal kernel. Additional numerical simulations show that the stability 
results remain valid not only for the infinitesimally small nonlinear modes, but also for the modes 
of finite amplitude. 
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I. INTRODUCTION 

Nonlocal nonlinearities are known to be of relevance 
for numerous physical applications. In optical appli¬ 
cations, nonlocal nonlinearity emerges in description of 
beam propagation in thermal media [I| or in atomic va¬ 
pors . Nonlocal nonlinearities have been shown to sup¬ 
port optical solitons M and to alter their properties, 
such as, for instance, stability Q, the collapse dynamics 
[3, or the interactions between solitons 

Nonlocality is also relevant in the theory of Bose- 
Einstein condensates (BECs) @ in view of the intrin¬ 
sically nonlocal nature of the inter-atomic interactions 
[S [TJ. Nonlocal nonlinearities are also inherent for 
plasma problems [T^ and for the theory on nematic crys¬ 
tals [T^ . 

In many relevant situations, nonlocal nonlinearity is 
accompanied by a certain potential which models a con¬ 
fining trap for a BEC or a refractive index modulation 
in an optical waveguide. Presence of the external poten¬ 
tial or nonlocality enriches the variety of possible non¬ 
linear modes allowing to observe the so-called higher 
(alias multi-pole) nonlinear modes 1,13 apart from the 
best studied fundamental (alias ground state) nonlinear 
modes [1, [3' Nonlocal nonlinearities have been consid¬ 
ered in combination with the double-well [l^ and the 
periodic potentials [iMl- It looks however that a sys¬ 
tematic study of nonlocal nonlinear modes in a parabolic 
potential has not been performed yet. In our present 
study, we make a step towards filling this gap and inves¬ 
tigate nonlocal nonlinear inodes in a complex parabolic 
T^T-symmetric potential of the form V{x) = — 2iax, 

where a is a real coefficient. For a = 0 the potential 
is reduced to the real parabolic potential whose modes 
are fairly well studied in the local case [3, [2^ , but have 
attracted much less attention in the nonlocal context. 

For a nonzero a, the potential V{x) is complex and VT 
symmetric [2l| , which is reflected by the fact that the real 


and imaginary parts of V(x) are even and odd functions 
of X, respectively. VT symmetry can be implemented 
in diverse physical setting, including optics [22|, atomic 
gases 3, and BECs [3|- 

Nonlinear modes and solitons have been studied in 
various PT-symmetric potentials with local [3]?7l| and 
nonlocal [3 nonlinearities, but the present work has a 
slightly peculiar focus. We are particularly interested 
in the examining the combined effect of the parabolic 
trapping, nonlocality and VT symmetry on the excited 
(higher-order) nonlinear modes. The focus of our study is 
additionally shifted towards the existence and stability of 
small-amplitude nonlinear modes, whose properties can 
be studied by means of asymptotic expansions describing 
bifurcations of nonlinear modes from the eigenstates of 
the underlying linear model. The linear spectrum of the 
VT -symmetric parabolic potential is available in the an¬ 
alytical, which makes it possible to reduce the stability 
problem to computing the real zeros of certain analytical 
expressions. Our analysis shows that either the degree of 
the nonlocality or the strength of the VT symmetry can 
be used to manage the stability of the excited modes. On 
the other hand, we reveal a dependence of the stability 
properties on the particular shape of the nonlocal kernel. 
We also touch upon the properties of nonlinear modes of 
finite amplitude and illustrate numerically that the sta¬ 
bility results hold not only for the infinitesimally small 
nonlinear modes but also remain valid in a finite range 
of the amplitudes of nonlinear modes. 

The paper is organized as follows. In the next section, 
we present a model to be considered throughout the pa¬ 
per. In Sec. imi we describe the linear properties of the 
model and introduce the families of small-amplitude non¬ 
linear modes. In Sec.|TV]we describe the methodology of 
the linear stability analysis. The main results of stabil¬ 
ity analysis are presented in Sec. |V] and Sec. |Vl] for the 
conservative and PT-symmetric cases, respectively. In 
Sec. IVIl] we demonstrate numerically that the presented 
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asymptotical stability results persist in a finite range of 
the amplitude of the modes. Concluding Sec. IVIIII out¬ 
lines the main results of the work. 


II. THE MODELS 


We use the dimensionless nonlocal nonlinear 
Schrddinger equation with an additional parabolic 
potential and nonlocal nonlinearity: 


.dq{x,z) d‘^q{x,z) ^ 


dz 


dx^ 

nco 


+ {x^ — 2iax)q(x, z) 


/ CO 

K{x — x')\q{x', z)\'^dx'. (1) 

-OO 


The case a = 0 corresponds to the real parabolic po¬ 
tential, while nonzero a governs the strength of the VT 
symmetry. The system is subject to gain in the domain 
X < 0 and is subject to losses for x > 0. In the op¬ 
tical context, equation m models beam guidance in a 
medium whose refractive index n(x) = nj.{x) + ini(x) 
has parabolic modulation of the real part nr{x) = x^, 
and linear modulation of the imaginary part rii = 2ax. 
Equation ([p with a = 0 also arises as a Gross-Pitaevskii 
equation [9| in the mean-field theory the harmonically 
trapped a quasi-one-dimensional BEC [13, in Hi]. In 
this case, q{x, z) is a quasi-one-dimensional wavefunction 
m , and z plays a role of time. Sign cr = ±I corresponds 
to attractive (focusing) and repulsive (defocusing) non- 
linearities, respectively. Nonlocal properties of the model 
are described by the kernel K (x) which is a positive func¬ 
tion with normalization f K(x)dx = 1 (hereafter all the 
integrals are taken over the whole real axis). Following 
to the previous studies on the subject H, H, [li, [ill , we 
pay the primary attention to the following two kernels: 

(i) The smooth Gaussian kernel 

K(x) = (2) 

VttA 

where A > 0 governs strength of nonlocality. 

(ii) Exponential kernel with the singularity at x = 0: 

(3) 


(iii) It will be also interesting to compare the ob¬ 
tained results with a slowly (algebraically) decaying ker¬ 
nel. While some studies (for example, [Ullil) deal with 
the so-called cut-off (CO) kernel which decays as |x|“^ 
as X —?► ±oo, in our study we choose the algebraically 
decaying Lorentzian kernel (the Lorentzian) 


K{x) 


X 


(5(x) 


if A > 0, 


if A = 0. 


(4) 


This choice is motivated by the fact that the Lorentzian 
kernel (jl]) allows for the relatively simple analytical eval¬ 
uation of the convolution integrals that arise in what fol¬ 
lows. 


III. BIFURCATIONS OF NONLINEAR MODES 

The stationary nonlinear modes for Eq. (HD admit the 
representation g(x, z) = e^^^w{x), where /3 is a real pa¬ 
rameter. It is convenient to introduce the shifted param¬ 
eter b defined through the relation /3 = b — which 
allows to rewrite the equation for the nonlinear modes in 
the following form: 

Wxx — bw — {x — ia)^w 

/ CO 

K{x — x')\w{x')\^dx'= 0, (5) 

-OO 

subject to the zero boundary conditions at the infinity: 
lima;^±oo w{x) = 0. 

Let us first recall the properties of the underlying linear 
problem which formally corresponds to cr = 0 in Eq. dSD- 
The resulting equation which is considered as a linear 
eigenvalue problem for b and w{x), has a purely real and 
discrete spectrum [13, [si, [3l| . The eigenvalues can be 

listed as = —(2n-|-l), n = 0,1,_ The corresponding 

eigenfunctions Wn{x) read 

Wnix) = CnHn{x - , (6) 

g-«V2 

Cn = , (7) 

v'V^2-n!G„(-2a2) 

where Hn(x) and Gn (x) are Hermite and Laguerre poly¬ 
nomials [ 33 , respectively. The choice of the coefficients 
Cn implies normalization J \wn\^dx = 1. Each Hermite 
polynomial Hnix) contains only powers of x of the same 
parity as the number n. This implies that in the con¬ 
servative limit (a = 0) the eigenfunctions are real-valued 
functions which have exactly n zeros and have the same 
parity as the number n. For o; ^ 0, each eigenfunction 
Wn{x) possesses real even part and odd imaginary part 
for even n and odd real part and even imaginary part for 
odd n. 

Returning now to the nonlinear problem, we look for 
bifurcations of small-amplitude nonlinear modes from the 
linear eigenstates. To this end, we introduce the formal 
asymptotic expansions which describe nonlinear modes 
Wn{x) branching off from the zero solution w{x) = 0 at 
the eigenvalues (27j |: 

ui„(x) = ewn+e'^w'^^'> +0{e^), K + +o{e^), 

( 8 ) 

where e 1 is a formal small real parameter (for the 
rigorous analysis of the bifurcations of small-amplitude 
modes see [^1. 

Since Wn{x) are normalized, the squared L^-norm of 
the nonlinear modes U = / |w„(x)pfix (which is also as¬ 
sociated with the total energy flow in the optical applica¬ 
tions of the model) is approximately equal to £^: 1/ ~ 
for e <C I. Substituting expansions ([5]) into Eq. and 
collecting the terms of the e^-order, one arrives at an in¬ 
homogeneous linear differential equation with respect to 
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(3) 

the function Wn ■ The solvability condition for the latter 

( 2 ) 

equation yields the expression for the coefficient ': 

W2) _ J wl{x) {J K{x - x')\wnix')\'^dx') dx 

JwUx)dx ■ 

For the conservative case {a = 0) the eigenfunctions 
Wn{x) are real, and hence the coefficients bn'^ are ob¬ 
viously real too for any n. For a ^ 0, the eigenfunctions 
Wn{x) possess nontrivial real and imaginary parts as de¬ 
scribed above. However, parities of real and imaginary 
parts of Wn (x) imply that for a real and even kernel func¬ 
tion K{x) = K(—x) = K*{x) the coefficients bn^ are also 
real, which shows that continuous families of nonlinear 
modes can exist in the PT-symmetric model with a > 0. 
The family with n = 0 represents the fundamental non¬ 
linear modes, and the families n = 1 , 2 ,... contain the 
higher-order modes. 

IV. LINEAR STABILITY OF NONLINEAR 
MODES 

A. Statement of the problem 

Following to the the standard procedure of the lin¬ 
ear stability analysis, we consider a perturbed nonlinear 
mode Wn{x) in the form q{x, z) = e’'^^[wn(x)+u(x)e^‘^^ + 
v*(x)e~'^‘^ which after linearization with respect to 
small functions u(x) and v(x) gives the linear stability 
eigenvalue problem 

Lp = ujp, p= , (10) 

with the linear operator L given as 

/ A crw„/Cn,„ 

\-o-w*/Cn,* -A^ 

where 


A = £ + a(bClWnl^) + (TWnbCw^, (11) 

= £1'-f cr(/C|w„p)-f (Tw*/C^„, (12) 

£ = d^ j dx'^ — b — {x — ia)"^, (13) 

£^ = d^/dx^ — b — {x + iaY, (14) 

K is the convolution operator, i.e., ICp = J K{x — 


x')p{x')dx', and JCp is the operator of “weighted” con¬ 
volution, i.e., K.pq = K{pq) = f K(x — x')p(x')q(x')dx'. 

The nonlinear mode Wn (x) is said to be linearly unsta¬ 
ble if operator L has at least one eigenvalue oj such that 
Im w < 0. Otherwise, Wn{x) is linearly stable. 

B. Linear stability of small-amplitude modes 

Let us now analyze the spectrum of the operator L for 
small-amplitude nonlinear modes belonging to the nth 


family and situated in the vicinity of the bifurcation from 
the linear limit. For nonlinear modes of zero amplitude 
(i.e., for e = 0) the operator L acquires the diagonal form 

^ = = -4)’ 

where 

(p 

£n =-^-bn - {x - (16) 

£i = ^-bn-(x + ia)^. (17) 

Since the spectrum of the operator £„ + bn is known 
from the discussion in Sec. imi it is easy to deduce that 
the spectrum of the operator £„ consists of two sub¬ 
sets. Eigenvalues and eigenvectors of the first subset 
read = 2(n - fc), = (wfe(a;),0)^, = 0,1,.... 

The second sequence reads = —2(n — k), p^^P = 
(0,w^(a;))'^, /c = 0,1,.... 

Thus, operator has a double zero eigenvalue Wnji = 
uJn,n = 0 which remains zero and double when pass¬ 
ing from the linear limit e = 0 to e > 0 [l^. The 
operator L„ also has 2n double nonzero eigenvalues: 
£ln,k = ‘^n'k = ^n^n-k^ where k runs from 0 to 2n except 
for k = n. Each double eigenvalue £ln,k is semi-simple; 
the two corresponding linearly independent eigenvectors 
are given as 

pS=('"o') PuSn-k=[^l_^)- (18) 

Passing from the linear limit e = 0 to small nonzero 
e, each double eigenvalue £ln,k generically splits into two 
simple eigenvalues which will be either both real or com¬ 
plex conjugated. If each double eigenvalue splits into two 
real ones, then the small-amplitude modes belonging to 
the nth family are stable for both attractive and repulsive 
nonlinearities. Vice versa, if at least one double eigen¬ 
value gives birth to a pair of two complex-conjugated 
simple eigenvalues, then the small-amplitude nonlinear 
modes of the nth family are unstable. Notice also that 
for n = 0 no double eigenvalues £tn,k exists. Therefore 
the fundamental solutions from family n = 0 are always 
stable in the small-amplitude limit. 

It is easy to check that the eigenvalue fln,k splits in 
exactly the same way as the opposite double eigenvalue 
£ln,2n-k = —^n,k-, i-G., if £ln,k Splits iiito two real eigen¬ 
values then Vln^ 2 n-k also splits into two real eigenvalues 
and vice versa. Therefore, in order to examine stability of 
small-amplitude nonlinear modes that belong to the nth 
family, it is sufficient to analyze only n positive double 
eigenvalues fln,k with A: = 0,..., n — I. 

In order to examine splitting of the double eigenval¬ 
ues £ln,k, we employ Eqs. ([5]) which yield the following 
asymptotic expansion for the linear stability operator: 
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L = Ln + ae^Ln'^ + o(e^), where the correction in ^ is 
given as 

(-hn'^ + (/Cllhnp) + WnJCw^ WnJC^^ \ 

y -w*K.^* - (/C|w;„p) - w*/CiD„/ 


Following the standard arguments of the perturbation 
theory for linear operators [s^, in order to address the 
splitting of the double eigenvalue iln,k we consider a 2 x 2 
matrix defined as 


Mn_k = 


P„,fc / 

/r(2)r,a) V 

\-On P„,fcJ Pn,2n-fe/ 

\ \^n,2n—k'> ■»^n,2n —fc/ 


(i 


D 

■t*n,2n—fc’ 


\ 


\Pn,fe’ Pn,A: / 


/r( 2 )r.(^^) n*'-'-'''* \ 

\^n Pn,2n-/c’Pn,2n-fc/ 

\Pn,2n—fe’ Pn,2n —fc/ / 


(//)* 


where (a, b) = J {x)a{x)dx for any two column vectors 
a and b (here ^ means the transposition and complex 
conjugation). If both the eigenvalues of matrix are 
real, then the double eigenvalue iln,k splits into two real 
and simple eigenvalues. If such a situation takes place for 
all fc = 0,1,..., n — I, then the small-amplitude nonlinear 
modes Wn(x) from the nth family are stable. On the 
other hand, if for some k the matrix Mn^k has a complex 
eigenvalue, then the double eigenvalue Vtn^k splits into a 
pair of complex-conjugated eigenvalues. This is sufficient 
to conclude that the small-amplitude nonlinear modes 
from the nth family are unstable. 

Recall that real and imaginary parts of Wn(x) have op¬ 
posite parities. Using this fact, one can establish that all 
the entries of the matrices Mn,k are real numbers which 
read 


{Mn^k)i_i — —b^n'^ + 

f Wk{u>kfCw*U>n + WnfCik^Wk} dx 

I wldx 

(Mn.k),^, = 6(2) - 

f xb2n—k{'dj2n—k^'w^'djn~\~'WnK^w^W2n—k^ dx 

f^L-kdx 


(19) 


( 20 ) 


(Adn,k)2^1 


f dx 

I^l^-kd-x 


( 21 ) 


{Mn^k)i2 


J '^n’^k^Wri‘^2n—k 

fwldx 


( 22 ) 


Thus in order to check reality of the eigenvalues of the 
matrix Mn^k, it is sufficient to consider the discriminant 
of its quadratic characteristic equation, i.e., the quantity 
Dn,k = [(M„,fc)^ ^-H(M„_fe)2_2]2-4detM„,fe. For Dn,k > 
0 the eigenvalues of Mn^k are real while a situation with 
Dn,k < 0 corresponds to a pair of complex conjugated 
eigenvalues. The situation Dn^k = 0 requires a more 
delicate analysis as in this case behavior of the double 


eigenvalue lln,k is determined by the next terms of the 
asymptotic expansions. 

Using analytical expressions for the eigenfuctions 
Wn(x), the entries of matrices Mn^k and hence the de¬ 
terminant Dn^k can be computed analytically with a 
computer algebra program (at least, for several relevant 
shapes of the kernel function K{x)). The resulting ana¬ 
lytical expressions are too bulky to be presented in the 
paper. However, for an interested reader to have an idea 
of how the expressions for the determinants Hn.fe read, 
in Appendix El we provide two expressions for Di q and 
£> 2,0 obtained for the Gaussian kernel and a = 0. Com¬ 
plexity of the expressions for Dn,k increases drastically 
if larger numbers n or nonzero values of a are involved. 
Nevertheless, these expressions are available analytically 
and can be used for finding the domains of positivity and 
negativity of the determinants Hn.fe which are considered 
as functions of A and a. 


V. RESULTS FOR THE REAL PARABOLIC 
POTENTIAL 


Let us start presentation of our results with an impor¬ 
tant particular case a = 0, when PT-symmetric contri¬ 
bution is absent, and the model is conservative. First, it 
is relevant to notice that for a = 0 the eigenvalue prob¬ 
lem (uni) admits an exact solution with eigenvalue w = 2 
and the respective eigenvector given as 


/ XWn + \ 

-XWr, +w'^ J 


(23) 


Notice however that this exact solution does not hold for 
nonzero a. 

Results of linear stability analysis of small-amplitude 
nonlinear modes from several first families for the Gaus¬ 
sian kernel m are summarized in Table HI For n = 
1, 2,... 6 and for k = 0,1,..., n — 1, Table U reports 
sign of the determinant Dn^k in the local limit A = 0 
and the approximated values of positive zeros of deter¬ 
minants Dn^k being considered as functions of A [at the 
same time, those are the values of A at which functions 
Dn,k = Dn,k{^) change their signs]. From Table [T] one 
can recover information on stability or instability of the 
small-amplitude nonlinear modes bifurcating from the 
nth linear eigenstate for the given value of the nonlocal¬ 
ity parameter A: the stability takes place if iAn,fe(A) > 0 
for each A: = 0,1,..., n — 1. 

It follows from Table U that for n = 1 and A: = 0 
the respective determinant iAi.o(A) is positive for any 
A > 0 (see also the explicit expression (EB) for Di o in 
Appendix . Therefore, the small-amplitude nonlinear 
modes bifurcating from the linear eigenstate n = 1 are 
stable for any A > 0. The small-amplitude modes with 
n = 2 are unstable for 0 < A < 1.82 but become stable for 
A > 1.82. For n > 2 the stability situation might be more 
complex with the half-axis A > 0 divided into several al¬ 
ternating intervals of stability and instability. In the local 
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TABLE I: Results of linear stability analysis for 0 = 0 and 
the Gaussian kernel. 

n k sign zeros of (if any) 


1 0 + 

2 0 

1 + 


no zeros 
1.82 

no zeros 


3 0 

1 

2 + 


2.36 

2.04 

no zeros 


4 0 + 

1 
2 

3 + 


2.61 2.85 

2.64 
2.26 

no zeros 


5 


0 

1 

2 

3 

4 


+ 


+ 


3.23 3.25 

3.12 
2.90 
2.45 

no zeros 


6 0 + 

1 + 

2 

3 

4 

5 + 


3.566 3.568 

3.45 3.51 

3.37 
3.13 
2.64 

no zeros 


case (A = 0) all the higher-order small-amplitude mode 
are unstable for n > 2 [ 2 ^. However, from Table |T] one 
observes that for each n there exists a threshold value 
An such that the small-amplitude modes become stable 
for the sufficiently strong nonlocality A > A„. The spe¬ 
cific values of A„ can be recovered from Table HI Ai = 0, 
A 2 Ri 1.82, A 3 Ri 2.36, etc. Thus the sufficiently strong 
nonlocality essentially alters (enhances) stability of the 
excited modes. 

We additionally mention that for each n one has 
.D„_„_i(A) > 0 for all A. This is a consequence of the 
exact solution (1231) which implies that the spectrum of op¬ 
erator L always contains an eigenvalue equal to 2. There¬ 
fore, the double eigenvalue = 2 always splits in 

the pair of real eigenvalues (one of which is equal to 2 ). 
As a result, the eigenvalue Q.n,n-i can not cause insta¬ 
bility of small-amplitude nonlinear modes. 

Let us now turn to examination of other kernel func¬ 
tions K{x). The results of linear stability analysis for the 
exponential kernel ([3|) are reported in Table HIl which fea¬ 
tures an essentially different stability picture comparing 
to that for the Gaussian kernel. Specifically, we observe 
that for the exponential kernel in the limit of strong non¬ 
locality A —>■ 00 , only the small-amplitude modes with 
n = 0,1,2 are stable, whereas the modes with n > 3 are 
unstable both in the local limit A = 0 and for all positive 


A. 

To summarize the intermediate results obtained so 
far, we have encountered different stability situations for 
two different shapes of the kernel K{x). In the case of 
the smooth Gaussian kernel, small-amplitude nonlinear 
modes bifurcating from the nth linear state become sta¬ 
ble in the limit of high nonlocality for any n. On contrary, 
for the exponential kernel the small-amplitude modes are 
stable in the strongly nonlocal limit only for n < 2 and 
unstable for n > 3. The found difference in stability 
situation of strongly nonlocal modes reveals certain cor¬ 
respondence to the previous results on the stability of the 
high-order nonlocal nonlinear states, such as multipole¬ 
mode ID solitons @ and 2D spatial solitons Both 
Q and 0 report that stable higher-order states can be 
found in the case of the Gaussian kernel, while the ex¬ 
ponential kernel results in instability of all high order 
solitons, starting from some threshold order. It has been 
also suggested in @ that the mentioned difference be¬ 
tween Gaussian and exponential kernels is related to the 
singularity which the exponential kernel has at the ori¬ 
gin a; = 0. Our approach allows to substantiate the above 
suggestion. Indeed, for the Gaussian kernel one can es¬ 
tablish asymptotic behavior of entries of the matrix Mn,k 
in the limit A —> 00 . Specihcally, for A —>■ 00 one has 

{Aln,k)l,l,{AIn,k)2,2, OC A“^, (24) 

(M„,fe) 2 . 1 ,(M„,fc)i ,2 cx a-( 2—2fc+i). (25) 

Technical details on computing the asymptotics (l24l) - 
(1^ are presented in Appendix [D1 Additionally, in Ap- 
pendices|B]and[Owe provide auxiliary relations necessary 
for the derivation of (IM 1 ) - (I^ . 

Relations ((Ml) - (l^5|) immediately imply that in the 
limit of large A for k = 0,l,...,n — 2 the diagonal el¬ 
ements of the matrix Mn^k dominate the off-diagonal el¬ 
ements, which guarantees reality of the eigenvalues of 
Mn,k- Since the instability never occurs for k = n — 1 
[due to the presence of the exact solution (l23l) . as dis¬ 
cussed above], we confirm that for the Gaussian kernel 
(or, more generally, for an infinitely differentiable kernel) 
the small-amplitude modes of arbitrarily high order n are 
stable in the limit of the strong nonlocality. 

At the same time, asymptotic results (I24l) - (l25l) . gener- 
ically speaking, do not hold for a non-differentiable ker¬ 
nel. For example, for the exponential kernel we have 
computed that for n = 3 and A: = 0 all the entries of 
the matrix Mn,k behave as oc for A —>■ oo, i.e., the 
diagonal domination in the strongly nonlocal limit is not 
guaranteed. 

Table uni which reports the stability analysis for the 
algebraically depending kernel (|T|) further confirms the 
stability of strongly nonlocal high-order small-amplitude 
modes for an infinitely differentiable kernel. Besides of 
this, we have also performed the same linear stability 
analysis for several nonphysical but still illustrative non¬ 
local kernels such as a rectangular kernel with finite sup- 
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TABLE II: Results of linear stability analysis for a = 0 and 
the exponential kernel. 

n k sign zeros of Dn.fe(A) 


1 

0 

-b 

no zeros 

2 

0 

- 

1.31 


1 

-b 

no zeros 

3 

0 

— 

no zeros 


1 

- 

1.54 


2 

-b 

no zeros 

4 

0 

-b 

7.25 


1 

- 

no zeros 


2 

- 

1.70 


3 

-b 

no zeros 

5 

0 

-b 

no zeros 


1 

- 

no zeros 


2 

- 

no zeros 


3 

- 

1.84 


4 

-b 

no zeros 

6 

0 

-b 

no zeros 


1 

-b 

no zeros 


2 

- 

no zeros 


3 

- 

no zeros 


4 

- 

1.97 


5 

-b 

no zeros 


port 

K{x)=l ^ fora:e(-A,A) ^^6) 

0 otherwise, 

a triangular shaped kernel with finite support 

f ^ + X for a: G (—A, 0) 

K{x) = <^ ^ + 1 for a; G (0, A) (27) 

I 0 otherwise, 

and a smooth exponentially decaying kernel 

X(X) = (ertc (^) + e,fc (i) 

(28) 

[the particular choice of the kernels (IMll - (|28)) is moti¬ 
vated by convenience of analytical evaluation of the re¬ 
spective convolution integrals]. Kernels (l26l) and (l28l) 
also feature stabilization of all higher-order strongly non¬ 
local modes, whereas for the triangular kernel (1^ the 
higher-order strongly nonlocal modes are unstable. This 
allows to conclude that the main feature responsible for 
the stability is not the very presence of the singularity of 
the kernel, but the location of the singularity: the rectan¬ 
gular kernel (E51) features two discontinuities at a: = ±A 


TABLE III: Results of linear stability analysis for a = 0 and 
the Lorentzian kernel. 


n 

k sign U„,fe(A)|;^^, 

j zeros of Un.fe(A) 

1 

0 

-b 

no zeros 

2 

0 

- 

1.52 


1 

-b 

no zeros 

3 

0 

— 

2.70 


1 

- 

1.66 


2 

-b 

no zeros 

4 

0 

-b 

2.77 3.446 


1 

- 

2.96 


2 

- 

1.77 


3 

-b 

no zeros 

5 

0 

-b 

3.88 3.991 


1 

- 

3.74 


2 

- 

3.21 


3 

- 

1.87 


4 

-b 

no zeros 

6 

0 

-b 

4.39 4.42 


1 

-b 

3.98 4.28 


2 

- 

4.00 


3 

- 

3.44 


4 

- 

1.97 


5 

-b 

no zeros 


but its behavior is nevertheless similar to that of the 
Gaussian and Lorentzian. The singularity at a: = 0 (the 
exponential and triangular kernels) does not allow for the 
stabilization of the higher-order small-amplitude modes. 


VI. RESULTS FOR THE PT-SYMMETRIC 
PARABOLIC POTENTIAL 

The linear stability analysis conducted above can be 
naturally extended on the case of nonzero VT -symmetric 
modulation a > 0. In this case, the stability results can 
be conveniently represented in the plane (a, A). Then the 
axis a = 0 recovers the results of the previous section, 
whereas the horizontal axis A = 0 corresponds to the 
case of local T^T-symmetric parabolic potential studied 
in [23 ■ The origin, i.e., the point A = a = 0 is “the 
most standard” case of the real potential with the local 
nonlinearity [ 2 ^. 

Similar to the previous section, we focus our attention 
on the Gaussian, exponential and Lorentzian kernels. Ac¬ 
cording to results of our analysis, for n = 1 the small- 
amplitude nonlinear modes are stable for any values of A 
and a from the considered range. For the families with 
n > 2, however, the stability situation becomes more 
complex and the (a, A)-diagrams turn out to be divided 
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into the domains of stability and instability. Such stabil¬ 
ity diagrams for the families with n = 2 and n = 3 are 
presented in Fig. [TJ 

Turning to the local limit (A = 0), one can observe that 
for sufficiently large a the small-amplitude modes become 
stable. Notice also that in the limit A 0 all the consid¬ 
ered kernels tend to the delta-function S(x) and therefore 
display the identical behavior. However, in accordance 
with results of the previous subsection, the considered 
kernels feature different behavior in the limit a = 0. As 
one can see in the panels of the first and the third rows, 
for the Gaussian and Loretzian kernels, in the limit a = 0 
large values of A lead to stabilization of the modes both 
for n = 2 and n = 3. For the exponential kernel (pan¬ 
els in the second row) the stabilization occurs for n = 2 
only, while for n = 3 the modes remain unstable for any 
arbitrarily large A. However, this situation changes dras¬ 
tically when a nonzero VT symmetry a > 0 is brought 
into consideration: say, already for a = 0.5 the nonlin¬ 
ear modes for the exponential kernel and n = 3 become 
stable for any A (from the considered range). Overall, we 
can conclude that the diagrams presented in Fig. [T] fea¬ 
ture nontrivial and fairly interesting structure. In par¬ 
ticular, one can observe that in some cases increase of a 
can lead either to stabilization or destabilization of the 
modes. The same remains true if one considers increase 
or decrease of A. The stability diagrams allow to con¬ 
jecture that for any given A the small-amplitude modes 
eventually become stable in the limit a —>■ oo. However, 
in order to verify this conjecture, an additional study 
must be undertaken. 

Notice also that for the case a > 0 allows for instability 
caused by the eigenstate with k = n — 1 (i.e., fc = 1 for 
n = 2 and fc = 2 for n = 3). No such instabilities is 
possible in the case a = 0 due to the exact solution (051) . 


VII. NONLINEAR MODES OF FINITE 
AMPLITUDE 




FIG. 1: Linear limit stability diagrams for n = 2 and n = 3 
in the plane (a. A). Darker domains correspond to instability. 
Labels “fc = 0”, “fc = 1”, etc. at the unstable domains show 
the number of the double eigenvalue Q.n^k responsible for the 
instability. 


The stability results of the previous sections apply to 
nonlinear modes of infinitesimally small amplitude. One 
can check persistence of the stability results by construct¬ 
ing numerically the families of nonlinear modes of finite 
amplitude and computing the eigenvalues of the corre¬ 
sponding linear stability problem. The families of non¬ 
linear modes can be visualized in the form of continuous 
curves on the plane aU vs. b where U = f |w„(x)pdx is 
the total energy flow. Then the small-amplitude modes 
are situated in the vicinity of the points = —(2n -|-1), 
aU = 0. 

Let us first recall the results related to the local case 
(A = 0). The corresponding diagrams for the conser¬ 
vative case [1^ and for the PT-symmetric case 271 are 
presented in the two upper panels of Fig. [2] The panels 
show several families that bifurcate from the lowest lin¬ 
ear eigenvalues n = 0,1,2,3. In the conservative case, 
only the two lowest families {n = 0,1) are stable in 


the small-amplitude limit, whereas the small-amplitude 
modes with n = 2,3 are unstable, In the PT-symmetric 
case, all the four shown families are stable in the linear 
limit, but at least three of the four families loose stability 
as the strength of nonlinearity becomes sufficiently large. 

The diagrams for the nonlocal case A = 4 (with the 
Gaussian kernel) are shown in the lower panels of Fig. [2] 
According to Table |T] and Fig. [U at the chosen value 
of the nonlocality parameter, the nonlinear modes with 
n = 0,1, 2, 3 are stable in the small-amplitude limit. The 
numerical results in Fig.[2]indicate that the stability per¬ 
sists at least for small and moderate strengths of the non¬ 
linearity {\crU\ < 40). Notice also that in the panel with 
(a = 1,A = 4) the nonlocality cancels the merging be¬ 
tween the subsequent families which can be observed in 
the panel with (a = I, A = 0). We however emphasize 
that an accurate description of the limit of strong non¬ 
linearity \(jU\ —>• oo should be elaborated by means of 
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VT symmetry results in a nontrivial stability picture for 
high-order small-amplitude nonlinear modes. The spe¬ 
cific results of our work can be outlined as follows: 

• A set of continuous families of nonlinear modes 
exists in the nonlocal PT-symmetric nonlinear 
Schrodinger equation. The spectrum of the cor¬ 
responding linear eigenvalue problem is available 
in the analytical form, and thus the stability prob¬ 
lem for the small-amplitude nonlinear modes can 
be reduced to the searching real roots of certain 
analytical expressions. 

• In the conservative case, small-amplitude nonlin¬ 
ear modes of arbitrary order become stable for suf¬ 
ficiently strong nonlocality, provided that the ker¬ 
nel is a sufficiently smooth function in the vicinity 
of the origin. If the kernel feature a singularity 
at a; = 0, no stabilization of high-order nonlinear 
modes is observed. 


FIG. 2: The families of nonlinear modes with n = 0,1,2, 3 
for different a and A (the Gaussian kernel). Stable modes 
correspond to the bold fragments of the curves. 


• Both the degree of nonlocality and the strength of 
the V'T symmetry can be used to manage stability 
of the small-amplitude modes. 


a proper asymptotic technique. While in the local case 
such studies have been already initiated [sj, |35| , analysis 
of the nonlocal VT -symmetric case in the limit of strong 
nonlinearity remains the open issue at the moment. 


• The above stability conclusions remain valid for 
nonlinear modes of finite amplitude. 


VIII. CONCLUSION 


Appendix A: Examples of expressions for Dn,k 


In this study, we have demonstrated that the com¬ 
bined effect of the parabolic potential, nonlocality and 


Let us present expressions Di q and Z? 2 ,o computed for 
the Gaussian kernel and a = 0: 


Difi 

D2,0 


1 {2X^ + lf 
4 ^(A2 + 2)®’ 

1 13968 Ai® - 28032 A^^ - 584 A^^ - 177136 A^^ - 54255 A® - 7564 A® -b 4306 A^ - 948 A^ - 207 
M TT (A2 -b2)^® 


(Al) 

(A2) 


Obviously, is positive for any A, while Z? 2 ,o is nega¬ 
tive for A = 0, but becomes positive for sufficiently large 
A. Using the standard numerical tools for searching the 
roots of polynomials, we find that D 2 ^o has two real ze¬ 
ros: A ~ ±1.82. The positive one is presented in the 
respective row {n = 2, k = 0) of Table ID 

ff(s) = T{f} = j f(x)V^^dx, (Bl) 

/(x) = T-^{g} = ^ j (B2) 


Appendix B: Definition and some properties of the 
Fourier transform 

We define the forward and backward Fourier trans¬ 
forms as 
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Then the Plancherel theorem says that 

J fi{x)f 2 {x)dx = -^ J gi{s)g 2 {s)ds, (B3) 

where 31 , 2 ( 5 ) = and from the convolution theo¬ 

rem one has 


T{h2}=HK}Hhi}, 


(B4) 


where h 2 {x) = f K(x — x')hi(x') dx'. 

Another well-known relation expresses the pth central 
moment of a function via the Fourier transform {p = 
0 , 1 ,...): 


f dP 

J xPf{x) dx = {-i)P ^5(s) 


(B5) 


s=0 


Appendix C: Some properties of the Hermite—Gauss 
eigenfunct ions 

Hereafter we assume the conservative case a = 0. 
The Hermite-Gauss eigenfunctions Wn{x) are defined by 
Eqs. ©, see [13 for definition of the Hermite and La- 
guerre polynomials. For a = 0 the modes Wn{x) are 
real-valued and satisfy the normalization 


/ CO 

w^ix) dx = 1. 

-00 


(Cl) 


The eigenfunctions are mutually orthogonal, i.e., 
/ WnWk dx = 0 for any n and k ^ n. Moreover, for 
any n and k < n: 

J x^Wnix) Wkix) dx = 0, p = 0,1,... ,n — k — 1,(C2) 


x” '"wn{x)wk{x) dx = (C3) 


In particular. 


j xw„{x)w„-i{x) dx = (C4) 


It is also useful to notice that 


J x^wlix) dx = n-f i. 


(C5) 


Appendix D: Asymptotics in the strongly nonlocal 
limit 


We assume the conservative case a = 0 and an in¬ 
finitely differentiable kernel written down as 


K{x) = l/r„ (f) 


(Dl) 


where function Kq{x) is even, real, and infinitely differ¬ 
entiable for any x. Then one can adopt the Taylor series 


“ /v(p)('n')TP 


(D2) 


P—0 


For the sake of brevity, in what follows let us use the 
notation 

Fn,k{s) = F{WnWk} = j Wn{x)wk{x)F^'^dx. (D3) 

Asymptotics for We start from Eq. , use 

Eqs. ()B3p . (jB4|l . ()C1|1 . and ()D2p . as well as reality of the 

(o') 

function Fnn{s), and evaluate bn as follows: 

b^^^ = {2n)-^ j F{K}Fl^{s) ds = 

J K{x)F-^{Fl^{s)} dx = 


P—0 




p=0 


AP+V! 


dsP 


s=0 


(in order to obtain the latter equality we also used (IBSIl l. 

From the definition of Eri_„(s) and from Eqs. m- 
()C5[1 . one can establish the following expansion: 


Fn,n{s) = I - i -f i -f 0 ( 5 ^) 


s^O- 


(D4) 


Therefore 


= (D5) 

Asymptotics for {Mn,k)i,i- Since the asymptotics 
for bn'^ is known from Eq. (ID5I) . it is convenient to con- 

(O') 

sider the quantity {Mn,k)i,i-\-bn which can be evaluated 
as 


^ (O') 

E {Fn,n{s)Fk,k{s) + |■F’n.fc(s)P)^^g ■ 

p=0 

Using (lD4l) and that 

Fn,k{s) = is j xwn{x)wk{x) dx + o{s)s^o, (D6) 

one can find the resulting asymptotics 

ATg (0) (n — k + 2(f xWnWk dx)^) 


dP 


{]\dn,k)l,l — 


2A3 

-l-0(A '^)a^+oo. 


The terms decaying as A ^ cancel each other. 
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For k = 0,1,... ,n — 2 one has J xWnWk dx = 0 and 
therefore 

W^n,fcjl,l —- ^^^3 -^ jA-)-+oo, 

while for A: = n — 1 one can use Eq. ()C4(1 and obtain 


nr 1 ^3)(« + i) , nn-4. 

— - ^ - hUtA , 


\—^“|“00 • 


Therefore several first terms of the series (lD7l) are zero, 
and the first nonzero term corresponds to p = 2{n — k) 
(provided that (0) ^ 0). Specifically, 

\(2n-k)\ 

{Mn,kh,i = TT^k 


2n-fc((„ _ ^)!)2y2n-2fe+l 


+0 


1 


\ 2n - 2k-\-2 


k\ 


A—^“1“00 


Asymptotics for {Mn,k) 2 ,i- Following to the same 
ideas, we can write down 


(^Mn,k^2,l — 




dp 


^ AP+ip! 


i-^r^{Fn,2n-,is)F:^kis)),^^ 


s=0 


(D7) 

From definition of functions F„_fe(s) and Fn, 2 n-k{s) it 
follows that in vicinity of s = 0 they behave as 


Finally we notice that proceeding in exactly the same 
way, one can prove analogous behavior of the other two 
entries Mn^k- Namely, {Mn,k) 2,2 decays as A“^, while 
(A7n,fe)i,2 decays as for A —>■ oo. 


Acknowledgments 


^—k 




F2n-k,k{s) = 


(n — fc)! V fc! 2 ” ^ 


^ s"-'' + o(s”-'^') 


s —^0 ? 


•n —fc 


{n — k)\\ n! 2” 
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